#ifndef __DIAMOND_SQUARE
#define __DIAMOND_SQUARE

#include <iostream>
#include <cstdlib>

using namespace std;

double randomvalue(int max_value) {
	return (rand()%(max_value));
}

double** generate_heightmap(int DATA_SIZE, double randomness_range) {
	// An initial seed value for the corners of the data
	double SEED = 0;

	double **data = new double*[DATA_SIZE];
	for (int i = 0; i < DATA_SIZE; ++i)
		data[i] = new double[DATA_SIZE];

	// Seed the data
	data[0][0] = SEED;
	data[0][DATA_SIZE-1] = SEED;
	data[DATA_SIZE-1][0] = SEED;
	data[DATA_SIZE-1][DATA_SIZE-1] = SEED;

	// Side length is length of a single square side
	// or distance of diagonal in diamond

	// For each iteration, we do one square step and one diamond step.
	for(int sideLength = DATA_SIZE-1;
	    // Side length must be >= 2 so we always have
	    // a new value (if it's 1 we overwrite existing values
	    // on the last iteration)
	    sideLength >= 2;
	    // Each iteration we are looking at smaller squares and
	    // diamonds, and we decrease the variation of the offset
	    sideLength /=2, randomness_range/= 2.0){
		// Half the length of the side of a square
		// or distance from diamond center to one corner
		// (just to make calcs below a little clearer)
		int halfSide = sideLength/2;

		//generate the new square values
		for(int x=0;x<DATA_SIZE-1;x+=sideLength){
			for(int y=0;y<DATA_SIZE-1;y+=sideLength){
				//x, y is upper left corner of square
				//calculate average of existing corners
				double avg = data[x][y] + //top left
					data[x+sideLength][y] +//top right
					data[x][y+sideLength] + //lower left
					data[x+sideLength][y+sideLength];//lower right
				avg /= 4.0;
				
				//center is average plus random offset
				data[x+halfSide][y+halfSide] = 
					//We calculate random value in range of 2h
					//and then subtract h so the end value is
					//in the range (-h, +h)
					avg + (randomvalue(2*randomness_range)) - randomness_range;
			}
		}
		
		//generate the diamond values
		//since the diamonds are staggered we only move x
		//by half side
		//NOTE: if the data shouldn't wrap then x < DATA_SIZE
		//to generate the far edge values
		for(int x=0;x<DATA_SIZE-1;x+=halfSide){
			//and y is x offset by half a side, but moved by
			//the full side length
			//NOTE: if the data shouldn't wrap then y < DATA_SIZE
			//to generate the far edge values
			for(int y=(x+halfSide)%sideLength;y<DATA_SIZE-1;y+=sideLength){
				//x, y is center of diamond
				//note we must use mod  and add DATA_SIZE for subtraction 
				//so that we can wrap around the array to find the corners
				double avg = 
					data[(x-halfSide+DATA_SIZE-1)%(DATA_SIZE-1)][y] + //left of center
					data[(x+halfSide)%(DATA_SIZE-1)][y] + //right of center
					data[x][(y+halfSide)%(DATA_SIZE-1)] + //below center
					data[x][(y-halfSide+DATA_SIZE-1)%(DATA_SIZE-1)]; //above center
				
				avg /= 4.0;
				
				//new value = average plus random offset
				//We calculate random value in range of 2h
				//and then subtract h so the end value is
				//in the range (-h, +h)
				avg = avg + (randomvalue(2*randomness_range)) - randomness_range;
				//update value for center of diamond
				data[x][y] = avg;
				
				//wrap values on the edges, remove
				//this and adjust loop condition above
				//for non-wrapping values.
				if(x == 0)  data[DATA_SIZE-1][y] = avg;
				if(y == 0)  data[x][DATA_SIZE-1] = avg;
			}
		}
	}
	
	for (int i=0; i<DATA_SIZE; i++) {
		for (int j=0; j<DATA_SIZE; j++) {
			data[i][j] /= 75.0;
		}
	}
	return data;
}

#endif
